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Abstract 

Useful approximation formulae for radiation impedance are given for the reflection 
coefficients of both infinitely flanged and unflanged rigid-walled cylindrical ducts. 
The expressions guarantee that simple but necessary physical and mathematical 
principles are met, like hermitian symmetry for the reflection coefficient (identical 
behaviour of positive and negative frequencies) and causality for the impulse re- 
sponse. A non causal but more accurate expression is also proposed that is suitable 
for frequency-domain applications. The formulae are obtained by analytical and 
numerical fitting to reference results from Levine & Schwinger for the unflanged 
case and extracted from the radiation impedance matrix given by Zorumski for the 
infinite flanged case. 
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1 Introduction and results 



The problem of the radiation acoustic impedance of the planar mode in a 
circular pipe with rigid walls (i.e. with homogeneous Neumann boundary con- 
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dition) is a classical problem of acoustics. Various detailed calculations for 
both unflanged and infinitely flanged ducts have already been provided, see 
e.g. [T|[2l5iri] . In addition, experimental investigations with various flanges have 
been compared to theoretical and numerical results [5]. There is a need for 
approximated formulae, such as those given by [lf6] or [TJ. Unfortunately these 
formulae do not fulfil the conditions for a physically representative model, as 
for instance the hermitian property of the reflection coefficient (see e.g. [8]): 

K(-u) = K*{oj) (1) 

where * means " complex conjugate of , or the causality of the impulse response 
of the reflection coefficient, obtained by inverse Fourier transform and denoted 
r(t) = FT -1 (lZ(u)), also known as the reflection function. This is essential to 
ensure that the time domain signals are real causal quantities. 

The aim of this short paper is to provide suitable approximate formulae based 
on accurate analytical or numerical fitting. Choosing the end of the bore as 
the reference plane and assuming a time dependence exp (—jut), the principal 
results will be expressed in terms of the pressure reflection coefficient 1Z(ou): 

n = -\n\e 2jkL = ^—±. (2) 

The dimensionless (i.e. divided by the characteristic impedance pc/Va 2 ) radi- 
ation impedance Z r can then be expressed in terms of TZ as follows: 

1-K J V 2 1 1 

k = ujjc and L denote the acoustic wave number and the end correction due 
to radiation, respectively. 

Section [2] describes how reference values were calculated. The requirements 
are discussed in Section |3j The approximate formulae are presented together 
with numerical results in Section HI Table [T] summarises the various formulae 
in the frequency and time domains. 



2 Calculation of the reference values 

2.1 Unflanged case 

Air vibrations in a rigid cylindrical pipe of negligible wall thickness and radius 
a propagate into free space through a circular sharp-edged opening. Nonlinear 
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effects like shock wave, vortex shedding, or mean axial flow are not consid- 
ered. Using the Wiener-Hopf technique, Levine & SchwingerJT] obtained an 
integral formulation of the reflection coefficient TZ(ka) for the planar mode, 
for frequencies below the cutoff frequency of the first non planar axisymmet- 
ric mode. We performed numerical evaluations of Eq. (V.16) of Ref. [JJ using 
the quadv function from Matlab[9]. The values obtained were compared to 
the asymptotic approximations given by Eqs. (VII. 1) and (VII. 2) of Ref. [JJ, 
with a maximum deviation compatible with the value of 3% mentioned in that 
paper. 



2.2 Flanged case 



The reflection coefficient for an infinitely flanged pipe has been calculated by 
Nomura et al [2] using Weber-Schafheitlin integrals and by Norris & Sheng[4j 
using modal expansion of the duct pressure field and a Green's function rep- 
resentation. The latter results in a modal sum for the planar mode reflection 
coefficient, each coefficient having a complex integral expression. 

It is also possible to derive the radiation impedance for the planar mode from 
the radiation impedance matrix written as in Ref. [10], considering only the 
modes being symmetrical about the axis of the pipe. The pressure and the 
velocity fields at the end of the bore are expressed in terms of the duct modes: 

P(r,k) = pc 2 Y,PnMr,k) and V(r, k) = c £ V n ^ n (r, k), (4) 

n>0 n>0 

with p the mean air density, and c the sound velocity in the free space, re- 
spectively. The dimensionless coefficients P n and V n are linked by P = ZV 
with: 

Z n , m = ^ / / /. fc ' r F V»(r, fc)V'm(r', k)dS{v)dS{v') Vn, m G [0, N], 

Z7T J Jr,r' Z7C\V — T \ 

(5) 

where the ip n { r ) are the normalised Bessel basis functions in the duct as in 
Eq. (8) of Ref. [TO]: 

ka \Jo{jn)\ 

] n being the nth zero of the Bessel function J x . Computation of the impedance 
matrix has been performed using Eq. (24) of Ref. [10J, giving results equal to 
the ones given by formulation given by Eq. (13) of Ref. [4]. 

The vectors P and V and the matrix Z are then decomposed into blocks, 
separating planar and non-planar evanescent components of the pressure and 
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Fig. 1. Real and imaginary parts of the radiation impedances of the circular piston 
set in an infinite plane baffle (zoo, dash-dotted curves) and of the planar mode of a 
cylinder in the flanged (Z r , solid curves) and unhanged cases (dashed lines). 

the flow fields: 

Po 
P' 

where zqo is the dimensionless radiation impedance of a circular piston and z 
gives the coefficients of the higher-order modes for an incident plane velocity 
field V(r) = Vq. Provided that the non-planar modes do not propagate, i.e. 
ka < ji ~ 3.832, and that the upstream bore is sufficiently long so that the 
evanescent modes do not meet a reflective obstacle in the duct, the following 
relation links the non-planar components: 

P' = -Z' C V (8) 

where Z' c is the diagonal matrix of the dimensionless characteristic impedances 
Z' c ,n = k/k n , with k n a = — (ka) 2 the wave number of evanescent modes. 
Eventually, planar components of pressure and flow are such that: 

Po = z 00 Vo + z T V = (^ o - z T (Z' + Z'J^z) V = Z r V (9) 

this expression showing the influence of evanescent non-planar modes present 
at the opening on the effective plane radiation impedance. The reflection coef- 
ficient for the planar mode is then obtained by means of Eq. ([2]). Fig. [1] shows 
the effects of the flange (dashed vs solid lines) and of the production of higher 
components at the end of the bore (solid vs dash-dotted lines). 
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The numerical method described here has been compared to the first coeffi- 
cient -Rooo °f the generalised reflection coefficient given by Zorumski (Eq. (32) 
of Ref. [TU]) and to the method given by Norris (Eqs.(ll) and (14) of Ref.Jl]) 
where the a nm coefficients can be deduced from the radiation impedance ma- 
trix. Results are in very good agreement, with relative error less than 10 -6 . 
Convergence of the calculations with respect to the number of higher order 
modes taken into account appeared to be significant for 10 modes. Matrix 
calculations were done using LAPACKjll], a standard linear algebra library. 



3 Model requirements 

Norris & Sheng (see Eq. (27) of Ref. [I]) attempted to provide simple accurate 
formulae for \1Z\ with a rational function. But a non negligible drawback of 
these formulae (and of other low-frequency approximations, see e.g. [6|7] ) is 
that they are not hermitian, i.e. the property TZ(— uj) = 1Z{uj)* is not satis- 
fied, and cannot be used for both positive and negative frequency domains. 
Another point is that it was not possible to estimate a time-domain reflec- 
tion function r(t) because the approximated formulae were not applicable for 
higher frequencies (the expression for the modulus |72-(o;)| become negative 
when uj increases). In this section the complete set of requirements on TZ(uj) 
is addressed. 



3. 1 Low frequency behaviour of the reflection coefficient 

The first requirement on the reflection coefficient is satisfaction of the following 
asymptotic forms in the low frequency domain: 



using the little-o notation. In the case of the baffled circular piston, the trun- 
cated part of the expansion is known to grow as (ka) 4 , but is more complicated 
in the unflanged case. Concerning the length correction as defined in Eq. ([2]), 
it becomes: 



K\(u 



0) = 1 - (3(kaf + o{ka) 2 with (5 





for the unflanged case 
for the flanged case, 



(11) 



these static values being the ones given by [T] and [I]. 
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3.2 Impulse response 



In order to enable the use of a time-domain reflection function, the condi- 
tion for inverse Fourier transform existence has to be fulfilled. First of all, 
|72.(k>)| has to tend to zero at high frequencies. The impulse response r(t) 
must also be a real quantity, which is guaranteed by the hermitian property 
in the frequency-domain (see e.g Eq. [[]). 

Another physics-driven requirement is that the impulse response should be 
causal. The Kramers- Kronig relations |12] provide a means to express this re- 
quirement in the frequency-domain. Another method is to study the placement 
of the poles of the reflection coefficient. Using the time dependence exp (—jut), 
the poles have to be located in the complex half-plane Im(cj) < (see e.g. [8], 
with the opposite time dependence convention); otherwise anticausal compo- 
nents appear in the inverse Fourier transform. However, it is not possible to 
calculate a time-domain impulse response based on the equations described 
in the previous section as this computation requires the knowledge of TZ(u>) 
on the full frequency range. As a practical alternative, we aim to provide ap- 
proximations having an extended domain of validity in both frequency and 
time. 

We notice that if the reflection coefficient is causal, the input impedance is 
causal as well, because the modulus of TZ(u>) is less than unity. Expanding 
(1 — TZ)^ 1 in Eq. ([3]) in the form of 1 + TZ + TZ 2 + . . ., the inverse Fourier 
transform of Z r is expressed as an infinite series of terms involving convolution 
products of causal functions, thus the result is a causal function. A physical 
interpretation in terms of successive reflections is classical. This reasoning can 
not be done in a reciprocal manner using Eq. (T5]), since the modulus of the 
dimensionless impedance can be large. 



4 Approximate formulae and results 



Two models are presented in this section. Both models satisfy the requirements 
mentioned in Section [3j they are causal in the time domain and satisfy the 
hermitian property and the desired low-frequency behaviour in the frequency 
domain. In particular they both result in very similar accuracy in the frequency 
domain below the below cutoff frequency ka = j x - The main difference is in 
the high-frequency behaviour which results in different smoothness properties 
at t = 0. 
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4.1 Model (u, a) 



The impulse response is modelled by the following expression: 

/ ' ct\ v ct 
r(t) = -A — exp (-a—) for t > 0, otherwise. (12) 
V a J a 

This definition ensures causality and hermitian symmetry (with real a and v). 
Its frequency domain representation (with the convention exp — jut) is given 
by: 

The adjustment of the low-frequency behaviour allows the determination of 
the parameters A, v and a: 

a ca u+1 ii , 2r? 2 

A = — r, a = -j- and i/ + 1 = — '— (14) 

or(i/ + i)' /3 /? v 7 



The value of A leads to the following expression : 

*M = - ( i - —Y ( " +1) ■ (is) 



a 



Numerical values of the parameters are given in Table [T] and frequency and 
time responses are displayed in Fig. [2] showing good agreement in frequency 
domain up to ka = 2 and monotonic decrease above. The time-domain reflec- 
tion function as defined by Eq. (fT2"j) is causal, as expected. 



4-2 Low order Pade approximant 



In order to permit analytical simplicity and possible calculations in the fre- 
quency domain, it is convenient to assume that the coefficient reflection de- 
pends on integers powers of frequency. That is not the case for the model (u,a) 
which has a non-integer exponent v. 

A rational function approximation of the form: 

. . 1 — riijka 

^ ~ ~l-d 1 jka + d 2 (jka) 2 ^ ' 

has been fitted to reference values. Low-frequency behaviour constrains the 
parameters nx, d\ and d 2 according to the following relations: 

d l -n l = 2rj, and d\-n\- 2d 2 = 2(3, (17) 
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Fig. 2. Comparison in frequency domain of the references (dashed lines) and model 
{v,a) of Eqs. (|12p through (|15p (in solid lines) for the unflanged (left column) 
and infinitely flanged (right column) cases. (a,b): modulus \TZ(ka)\ of the reflection 
coefficient (approximation errors in solid lines in (e,f)). (c,d): dimensionless length 
correction L/a (approximation errors in dash-dotted lines in (e,f)). (g,h): example 
of the reflection function (a = 7mm). 



so that a one-dimensional fitting is done using a Nelder-Mead simplex algo- 
rithm to minimise the error between the complex evaluations of the reference 
and the model. N = 75 values of TZ(ka) were computed for regularly spaced 
values of ka between and ji (the cutoff value for the first axisymmetric 
non planar mode). Numerical values obtained by constrained optimisation 
(requiring causality or, equivalently, that no pole lie in Im(u;) > 0) are given 
in Table[TJ In the same manner as for the previous model, frequency and time 
responses are shown in Fig. [3j TZ(uj) seems almost as well adjusted to the ref- 
erence values as the model (z/,a) but, unlike that model, the impulse response 
r(t) now shows an instantaneous initial step at t — 0: 

r ( f > °) = ~TT~ V - l)e" 72 ^ - (n l7 i - l)e~^f) , (18) 

ad 2 {jx - 72) v 1 

where 

7i,2 = ^ (dx ± \jd\ - 4rf 2 ) • (19) 
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Fig. 3. Comparison in frequency domain of the references (dashed lines) and model 
(1,2) of Eqs. (|16p through (|20p (in solid lines) for the unhanged (left column) and 
infinitely flanged (right column) cases. (a,b): modulus \TZ(ka)\ of the reflection co- 
efficient (approximation errors in solid lines in (e,f)). (c,d): dimensionless length 
correction L/a (approximation errors in dash-dotted lines in (e,f)). (g,h): example 
of the reflection function (a = 7mm). 

The radiation impedance Z r is given by the following expression: 

(d 1 -n 1 )jka-d 2 (jka) 2 
r[U) ~2-(d 1 + ni )jka + d 2 (jkaf [ ' 

similar to the one suggested by Doutaut & Chaigne[13j. where the coefficients 
were obtained by numerical fitting on the approximated radiation impedance 
given by Causse et al.^\. 



5 Relaxing the causality constraint 



The formulae given previously may appear somewhat inaccurate especially as 
the cutoff frequency of the first higher order mode is approached. However, it 
should be kept in mind that they are intended to simultaneously approximate 
the low frequency behaviour of the reflection coefficient, satisfy the hermitian 
property, and produce a causal physical response. With all these constraints 
being satisfied, the relative errors on the modulus and the length correction 
are less than 8% for ka < 2. For applications where the time-domain response 
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Fig. 4. Comparison in frequency domain of the references (dashed lines) and the non 
causal model of Eqs. (|2ip and (|22p (in solid lines) for the unhanged (left column) 
and infinitely flanged (right column) cases. (a,b): modulus \1Z(ka)\ of the reflection 
coefficient (approximation errors in solid lines in (e,f)). (c,d): dimensionless length 
correction L/a (approximation errors in dash-dotted lines in (e,f)). 



is not a critical criteria, relaxation of the causality constraint may lead to 
improved approximations. 

We find that modelling \1Z{ka)\ and L{ka)/a with Pade approximants of order 
(2, 6) provides the more accurate approximations with the following expres- 
sions: 



\n\ 



1 + ai(ka) 2 
1 + (J3 + ai)(ka) 2 + a 2 {ka) i + a 3 (ka) 6 ' 
L_ l + b^kaf 

a~ V l + h 2 (ka) 2 + b 3 (ka) 4 + b A {kaf ' 



(21) 
(22) 



and the numerical values given in Table [TJ These formulae fulfil all of the 
requirements except the causality principle, and approximation errors for the 
modulus and the length correction are less than 2% for ka < 3 for both 
unflanged and infinitely flanged cases as shown in Fig. HI Combining formu- 
lae (12T!) and (1221) into Eq. ([3]) allows the computation of the acoustic radiation 
impedance. 
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6 Conclusion 



Practical approximation formulae of the acoustic radiation impedance of tubes 
are required for a variety of applications in acoustics. The necessity for such 
useful formulae is illustrated, in part, by the list of references cited here, which 
is long but certainly not exhaustive. The present results attempt to satisfy all 
the criteria required of the realistic physical system, such as causality, hermi- 
tian response and faithful low frequency behaviour. We therefore believe that 
the present contribution can be useful for calculations in the both time and 
frequency domains. Practical applications include, for instance, measurements 
in duct acoustics, musical acoustics, or loudspeakers enclosures where mean 
flow is absent or very slow. 

Difficulty remains for the digital (sampled time) domain. In fact, the approxi- 
mate formulae provided in this paper are all expressed in terms of the variable 
ha. A digital filter modelling the radiation at the end of a cylindrical duct has 
to be written as a function of the z = expjuT e variable (in the T e -sampled 
time domain) and the coefficient of the filter are then dependent on the ra- 
dius a. This could present difficulties if one wants to obtain approximation 
formulae suitable for the digital domain. 

Despite such remaining technical issues, the calibrated formulae summarised 
in Table [1] offer the acoustician practical and useful formulae for realistic ap- 
plications. The numerical results of the calculations described in Section [2] are 
available at the following url: |http : / / www . lma . cnrs-mrs . f r/~ PIM/DATA/ , 



Acknowledgements 

The study presented in this paper was conducted with the support of the 
French National Research Agency anr within the Consonnes project. The 
authors wish to thank Ch. Morfey and P. O. Mattei for useful discussion. 
A.N.N, wishes to express thanks to the Laboratoire de Mecanique Physique, 
Universite Bordeaux 1, for hosting him. 



References 

[1] H. Levine, J. Schwinger, On the radiation of sound from an unflanged circular 
pipe, Physical Review 73 (4) (1948) 383-406. 

[2] Y. Nomura, I. Yamamura, S. Inawashiro, On the acoustic radiation from a 



11 



Model 


Unflanged case 


Flanged case 




— I/O r> — f) 61 33 


B — 1 77 — O 8216 


KM = - (1 - i^y {l/+1) 

r(t) = -A(f) !y exp(-af) 


a = 1.226o 
jy = 0.504 
A — 1 ^/L c 


a = 0.8216 

A — n SRI c 




i 167 


1 i«2 


'M"- 1 ; l-d 1 jka+d 2 (jka) 2 








d 2 = 0.457 


d 2 = 0.649 




ai = 800 


ai = 730 




a2 = 0.26o 


02 = 0.3/2 


\-rpi l+ai(ka) 2 

~~ l+(/3+ai)(/ca) 2 +a 2 (fca) 4 +a 3 (fca) 6 ' 
L l+fei(fca) 2 
a — '/l+fe 2 (M 2 +MM 4 +MM 6 ' 


03 = 0.0263 
b\ = 0.0599 


03 = 0.0231 
Oi = 0.244 


non causal 


b 2 = 0.238 


62 = 0.723 




b 3 = -0.0153 


63 = -0.0198 




6 4 = 0.00150 


64 = 0.00366 


able 1 



This table summarises the results for the approximate radiation models correspond- 
ing to time dependency exp(— jut). For the opposite convention exp(+jcjt), com- 
plex conjugation of the expressions is needed. 

flanged circular pipe, Journal of the Physical Society of Japan 15 (1960) 510- 
517. 

[3] Y. Ando, On sound radiation from semi-infinite circular pipe of certain wall 
thickness, Acustica 22 (4) (1969) 219-225. 

[4] A. N. Norris, I. C. Sheng, Acoustic radiation from a circular pipe with an infinite 
flange, Journal of Sound and Vibration 135 (1) (1989) 85-93. 

[5] J.-P. Dalmont, C. J. Nederveen, N. Joly, Radiation impedance of tubes with 
different flanges : numerical and experimental investigations, Journal of Sound 
and Vibration 244 (3) (2001) 505-534. 

[6] P. O. A. L. Davies, J. L. Bento Coelho, M. Bhattacharya, Reflection coefficients 
for an unflanged pipe with flow, Journal of Sound and Vibration 72 (4) (1980) 
543-546. 

[7] R. Causse, J. Kergomard, X. Lurton, Input impedance of brass musical 
instruments — comparison between experiment and numerical models, Journal 
of the Acoustical Society of America 75 (1) (1984) 241-254. 

[8] R. W. Rienstra, Impedance models in time domain, including the extended 
helmholtz resonator model, in: 12th AIAA/CEAS Aeroacoustics Conference, 



12 



Vol. 2006-2686, Cambridge, MA, USA, 2006. 



[9] Matlab®, The MathWorks, Inc. 

[10] W. E. Zorumski, Generalized radiation impedances and reflection coefficients of 
circular and annular ducts, Journal of the Acoustical Society of America 54 (6) 
(1973) 1667-1673. 

[11] Lapack Users'Guide, Society for Industrial and Applied Mathematics, ISBN 
0-89871-447-8. 

[12] J. G. McDaniel, Applications of the causality condition to one-dimensional 
acoustic reflection problems, Journal of the Acoustical Society of America 
105 (5) (1999) 2710-2716. 

[13] V. Doutaut, D. Matignon, A. Chaigne, Numerical simulations of xylophones. 
II. Time-domain modeling of the resonator and of the radiated sound pressure, 
Journal of the Acoustical Society of America 104 (3) (1998) 1633-1647. 



13 



